Initialization

The urban heat island equations are based on Theeuwes (2017).

Initialization

wd<-"D:/uhimax/"

grd.svf<-stack(paste0(wd,"Grids_veg_svf/svf_utrecht_1m.grd"))
grd.fveg<-stack(paste0(wd,"Grids_veg_svf/greenness_utrecht_smooth_500m.grd"))
grd.fveg<-projectRaster(grd.fveg,crs=crs(grd.svf))

IUtrecht23<-data.frame(lon=52.079,lat=5.139)
coordinates(IUtrecht23)<- ~lat+lon
crs(IUtrecht23)<-crs("+init=epsg:4326")

mapview(grd.fveg) + mapview(IUtrecht23)
## Warning in rasterCheckSize(x, maxpixels = maxpixels): maximum number of pixels for Raster* viewing is 5e+05 ; 
## the supplied Raster* has 1348088 
##  ... decreasing Raster* resolution to 5e+05 pixels
##  to view full resolution set 'maxpixels =  1348088 '
IUtrecht23_fveg<-raster::extract(grd.fveg,IUtrecht23)
## Warning in .local(x, y, ...): Transforming SpatialPoints to the CRS of the
## Raster
IUtrecht23_svf<-raster::extract(grd.svf,IUtrecht23)
## Warning in .local(x, y, ...): Transforming SpatialPoints to the CRS of the
## Raster

From 10min data to Daily UHImax parameters

From the autmaomatic weather stations (AWS) in the rural area the following data is required:

The 10 minute data from global radiation, wind speed, temperature and pressure is used to calculate the UHImax. The hourly relative humidity, precipitation and mean wind are used to filter out synoptic situations with frontal system and fog. Under these different conditions there is no clear relation between the AWS and city temperature.

Cabauw_data<-prepare_Cabauw(ten_min = paste0(wd,"Cabauw_meteo/cabauw10min.csv"),
                            hourly_data = paste0(wd,"Cabauw_meteo/cabauw_hourly.csv"))
## Reading the 10min data
## Preparing 10min variables
## Calculating meteo parameters for Theeuwes(2017) using the 10min data
## Reading hourly data

Wunderground network data in Utrecht

Within the city of Utrecht we have carefully selected three stations from the Wunderground network:

These stations have been measuring the city temperature for the past couple of years. The data can be downloaded using download_time_seq, note that the temperature is in Farenheid and not degrees Celcius. The data has been filtered and interpolated to 10minute timestamps. Unrealistic high and low temperatures are set to NA. In case more than 8 identical measurements in a row occur the value is also set to NA.

Calculate UHI max parameters

Setting all parameters to required to calculate the UHI with the formula of Theeuwes (2017). For this figure no subset of meteorological conditions was made.

uhi_Utrecht<-uhimax_params(STN="UTRECHT23",
              svf=IUtrecht23_svf,
              fveg=IUtrecht23_fveg,
              city_T=paste0(wd,"Wunderground/Filtered/IUTRECHT23_filtered.rds"),
              rural_T=Cabauw_data$Cabauw_10min,
              rural_meteo=Cabauw_data$Cabauw_Theeuwes)
## Reading city temperature data
## Merging with rural temperature data
## Calculating the maximum UHI from temperature observations
ggplot(uhi_Utrecht,aes(Tref,Tcity,colour=S_new)) + geom_point() + geom_abline() + scale_color_gradientn(colours = heat.colors(20))

#Similarly the same function can be applied to other stations within Utrecht:
svf_grn<-readRDS(paste0(wd,"/Rdata/wunderground_svf_grn.rds"))
df.svf_grn<-data.frame(svf_grn)
df.svf_grn<-df.svf_grn[df.svf_grn$Station.ID %in% c("IUTRECHT196","IUTRECHT376","IUTRECHT299"),]
stations<-df.svf_grn$Station.ID
city_T<-paste0(wd,"Wunderground/Filtered/",df.svf_grn$Station.ID,"_filtered.rds")
svf<-df.svf_grn$svf
fveg<-df.svf_grn$grn

# uhimax_params(STN=stations[3],svf=svf[3],fveg=fveg[3],city_T = city_T[3],rural_T = Cabauw_data$Cabauw_10min,rural_meteo = Cabauw_data$Cabauw_Theeuwes)

Wunderground_Utrecht<-mapply(uhimax_params, STN = stations, svf = svf, fveg = fveg, city_T = city_T,
                      MoreArgs = list(
                           rural_T = Cabauw_data$Cabauw_10min,
                           rural_meteo = Cabauw_data$Cabauw_Theeuwes
                           ,wd = paste0(wd,"UHImax/"), write.file = TRUE
                           ),SIMPLIFY = FALSE)
## Reading city temperature data
## Merging with rural temperature data
## Calculating the maximum UHI from temperature observations
## writing file
## Reading city temperature data
## Merging with rural temperature data
## Calculating the maximum UHI from temperature observations
## writing file
## Reading city temperature data
## Merging with rural temperature data
## Calculating the maximum UHI from temperature observations
## writing file
Wunderground_Utrecht<-do.call("rbind",Wunderground_Utrecht)

ggplot(Wunderground_Utrecht,aes(Tref,Tcity,colour=S_new)) + geom_point() + geom_abline() + scale_color_gradientn(colours = heat.colors(20))

Subset based on weather conditions

uhimax_files<-list.files(paste0(wd,"UHImax/"),
                         full.names = TRUE,pattern=".txt")
uhimax_files<-lapply(uhimax_files,fread)
uhimax_Utrecht<-do.call("rbind",uhimax_files)

# uhimax_Utrecht<-fread(paste0(wd,"UHImax/IUTRECHT23_UHIparams.txt"))

uhimax_Utrecht<-merge(Cabauw_data$Cabauw_sub,uhimax_Utrecht,by=c("start","stop"))
uhimax_Utrecht<-uhimax_Utrecht[which(uhimax_Utrecht$Rain==TRUE),]
uhimax_Utrecht<-uhimax_Utrecht[which(uhimax_Utrecht$Wind==TRUE),]
uhimax_Utrecht<-uhimax_Utrecht[which(uhimax_Utrecht$rh==TRUE),]
# uhimax_Utrecht<-uhimax_Utrecht[which(uhimax_Utrecht$Select==TRUE),]
uhimax_Utrecht<-uhimax_Utrecht[which(uhimax_Utrecht$Tref>17),]
# uhimax_Utrecht<-uhimax_Utrecht[which(uhimax_Utrecht$stn %in% c("IUTRECHT23","IUTRECHT376")),]
uhimax_Utrecht<-uhimax_Utrecht[which(uhimax_Utrecht$stn %in% c("IUTRECHT23","IUTRECHT196","IUTRECHT376","IUTRECHT299")),]

Correlations with the three meteo params

The maximum urban heat island is in the cases we selected related to:

\[ UHI_{max} = (2-SVF-F_{veg}) \sqrt[4]{\frac{S*DTR^3}{U}} \] The station IUTRECHT23 has the best correlation coefficient with 0.6, while also including the other observations in the center of Utrecht reduces R2 to 0.4.

uhi_u<-ggplot(uhimax_Utrecht,aes(U,UHImeasured,colour=DTR)) +geom_point() + scale_color_gradientn(colours = terrain.colors(20))
uhi_DTR<-ggplot(uhimax_Utrecht,aes(DTR,UHImeasured,colour=S_new)) + geom_point()+scale_color_gradientn(colours = heat.colors(20))
uhi_S<-ggplot(uhimax_Utrecht,aes(S_new,UHImeasured,colour=U)) +geom_point() + geom_point()+scale_color_gradientn(colours = topo.colors(20))
ggsave(uhi_u,filename = paste0(wd,"UHImax/fig/uhi_u.png"))
## Saving 7 x 5 in image
ggsave(uhi_DTR,filename = paste0(wd,"UHImax/fig/uhi_DTR.png"))
## Saving 7 x 5 in image
ggsave(uhi_S,filename = paste0(wd,"UHImax/fig/uhi_S.png"))
## Saving 7 x 5 in image
uhi_u

uhi_DTR

uhi_S

ggplot(uhimax_Utrecht,aes(UHImeasured,(2-svf-fveg)*meteo,colour=factor(stn))) +geom_point() +geom_abline() +xlim(0,10)+ylim(0,10)
## Warning: Removed 1 rows containing missing values (geom_point).

caret::R2((2-uhimax_Utrecht$svf-uhimax_Utrecht$fveg)*uhimax_Utrecht$meteo,uhimax_Utrecht$UHImeasured)
## [1] 0.3996237
caret::RMSE((2-uhimax_Utrecht$svf-uhimax_Utrecht$fveg)*uhimax_Utrecht$meteo,uhimax_Utrecht$UHImeasured)
## [1] 1.095273

Using the formula of Theeuwes for Utrecht

The formula of Theeuwes(2017) was derived for:

Outside of this range the formula has not been tested.

svf.r<-raster::resample(grd.svf,grd.fveg,method="bilinear")
values(svf.r)[values(svf.r)<0.2 | values(svf.r)>0.95] = NA
values(grd.fveg)[values(grd.fveg)<0 | values(grd.fveg)>0.4] = NA

city_center<-extent(c(5.1,5.13,52.08,52.098))
city_center<-raster(city_center)
values(city_center)<-1
crs(city_center)<-CRS("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")
city_center<-projectRaster(city_center,crs=crs(grd.fveg))

st<-(2-svf.r-grd.fveg)*mean(uhimax_Utrecht$meteo)
st<-crop(st,city_center)
mapview(st)